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ABSTRACT 


This  thesis  examines  the  performance  of  three  different  orbital  propagators  to  determine 
which  provide  the  best  performance  for  use  in  Low  Earth  Orbit  Rendezvous.  The  performance 
evaluation  is  based  upon  the  propagator’s  accuracy  and  the  amount  of  time  required  to  produce  a 
solution.  A  Cowell  high-fidelity  propagator  is  used  as  a  base  line  for  comparison  with  an  Encke 
and  Clohessy-Wiltshire  propagator.  To  further  enhance  the  examination  a  Jacchia-70 
atmospheric  model  and  a  GEM-9  Geopotential  model  are  used  to  provide  perturbing  acceleration 
inputs  to  the  propagators.  All  comparisons  are  performed  in  a  Local  Vertical,  Local  Horizontal 
Reference  Frame  with  the  target  spacecraft  at  the  coordinate  center.  Tainting  of  the  input  data  set 
by  a  prior  processor  make  the  findings  suspect.  Findings  support  the  prediction  that  while  the 
Cowell  propagator  is  the  most  accurate  it  also  takes  the  most  time  to  achieve  results.  Also,  the 
Clohessy-Wiltshire,  while  taking  the  least  time  is  the  most  inaccurate.  The  Encke  propagator 
deliveries  the  most  balanced  result. 
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I.  INTRODUCTION 


The  intercept  of  two  large  objects  in  a  gravity  free  environment  can  either  be  very 
easy  or  very  hard  depending  on  whether  or  not  you  want  both  objects  intact  after  the 
meeting.  If  not,  a  high  energy  collision  works  every  time.  For  rendezvous  in  space,  the 
object  is  for  two  spacecraft  to  meet  in  space,  work  together  and  separate  with  no  damage 
to  either.  This  requires  a  great  deal  of  knowledge  and  planning  involving  a  medium  in 
which  human  beings  have  little  intuition.  To  counter  this  shortcoming  tools  have 
developed  to  help. 

Orbital  propagators  have  been  used  for  years  to  predict  the  positions  and 
velocities  of  objects  in  orbit.  They  are  an  integral  part  of  space  navigation.  Their  only 
real  shortcoming  is  their  size  and  complexity.  A  good  high  fidelity  propagator  can 
involve  thousands  of  lines  of  code  and  many  minutes  of  processing  time  in  order  to 
predict  a  spacecraft’s  orbital  location.  While  this  is  acceptable  for  mission  planning,  it  is 
not  very  useful  for  real-time  propagation. 

In  September  of  1993  NASA  launched  STS-51  which  experimented  with  a  new 
navigational  tool.  The  Detailed  Technical  Objective  700,  consisted  of  a  GRID  laptop 
computer  with  software  written  by  NASA  and  a  group  of  students  from  the  Naval 
Postgraduate  School  (Ref.  1).  As  a  precursor  to  an  automated  navigational  system,  this 
experiment  was  designed  to  demonstrate  how  on-orbit  orbital  prediction  and  graphical 
displays  could  assist  astronauts  navigating  of  the  Shuttle  Orbiter.  Such  an  aid  was  found 
to  be  especially  valuable  during  Proximity  Operations  and  Rendezvous  (Ref.  1).  The 
programs  RPOP  and  NPSLite  flew  and  were  used  during  an  actual  rendezvous. 

The  problem  with  RPOP  or  NPSLite  as  bases  for  an  automated  navigation 
program  lays  in  their  design.  RPOP  used  a  set  of  differential  equations  developed  by 
Clohessy-Wiltshire  for  orbital  prediction  while  NPSLite  used  a  “f  ’  and  “g”  function  type 
propagator  (Ref.  1).  Neither  of  the  propagators  took  into  account  perturbations  caused  by 
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external  forces  except  those  inherent  in  the  Two  Body  Problem.  This  was  not  an  arbitrary 
decision.  Normally  the  equations  used  to  determine  orbital  position  for  a  two  body 
problem  are  not  that  complex,  and  therefore  fairly  simple  for  a  computer  to  handle.  The 
addition  of  various  gravitation  and  atmospheric  models  to  more  accurately  reflect  small 
perturbing  accelerations  can  greatly  increase  both  program  size  and  processing  time.  In 
the  case  of  RPOP  and  NPSLite,  it  was  assumed  that  two  objects  in  nearly  identical  orbits 
are  perturbed  relatively  the  same.  Therefore,  when  dealing  with  a  purely  relative 
reference  frame  these  perturbations  could  be  ignored  with  little  loss  of  accuracy.  The 
purpose  of  this  thesis  is  to  determine  the  validity  of  that  assumption. 

As  a  result  of  the  success  of  RPOP  and  NPSLite,  NASA  has  decided  to  develop  a 
flight  qualified  version  of  the  RPOP  software.  It  is  NASA’s  intention  to  have  the 
software  available  for  use  during  the  future  Space  Shuttle  -  MIR  mission  in  1995.  One  of 
the  difficulties  in  the  flight  certification  is  verifying  that  the  algorithm  upon  which  the 
software  is  based  is  sufficiently  accurate.  By  expanding  this  query,  we  can  quickly  find 
ourselves  asking  which  of  the  many  orbital  propagation  models  and  algorithms  available 
is  the  best  for  the  purpose  at  hand.  This  is  the  second  objective  of  this  thesis. 

As  previously  mentioned,  the  on-orbit  computation  time  a  propagator  requires  can 
be  critical.  In  order  to  be  useful  for  the  flight  crew,  a  propagator  should  accurately  predict 
the  Orbiter’s  position  and  velocity  at  least  one  orbit,  approximately  90  minutes,  in 
advance.  Because  even  the  most  accurate  propagator  is  vulnerable  to  errors,  this  90 
minute  requirement  must  be  broken  up  into  a  series  of  steps  typically  on  the  order  of  a 
few  seconds.  This  in  turn  means  that  the  propagator  must  provide  thousands  of 
predictions  for  each  90  minute  propagation.  Additionally,  to  be  of  any  use  on-orbit,  the 
software  needs  to  perform  all  of  these  predictions  repeatably  every  few  seconds.  For  a 
high  fidelity,  time  intensive  propagator,  this  time  restriction  may  render  the  algorithm  of 
little  value. 
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Because  of  the  time  constraint,  the  best  propagator  is  not  necessarily  the  most 
useful.  The  best  propagator  balances  accuracy  with  low  computation  time.  An  ideal 
propagator  would  be  best  in  both  categories,  but  as  will  be  seen  this  is  rarely  if  ever  the 
case.  In  fact,  the  ideal  propagator  may  eventually  depend  on  far  more  than  just  accuracy 
and  speed,  but  for  the  moment  this  study  will  be  limited  to  just  these  two  criteria.  It  is  the 
intent  of  this  thesis  to  provide  the  initial  look  and  possibly  a  guideline  to  future  research 
in  the  ultimate  choice  of  a  best  orbit  propagation  model  and  algorithm. 
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II.  REFERENCE  FRAMES 


Which  reference  frame  a  body  is  being  observed  in  is  integral  to  defining  and 
describing  the  forces  interacting  with  it.  By  definition,  the  Newton’s  laws  of  motion 
apply  only  in  an  inertial  reference  frame.  In  order  to  translate  to  a  non-inertial  reference 
frame  and  use  the  same  laws,  one  must  translate  those  laws  into  an  appropriate  form. 

Integral  to  the  definition  of  a  reference  frame  is  the  coordinate  system  used  to 
define  location  within  said  frame.  Basically,  all  a  coordinate  system  does  is  determine  the 
location  of  a  point  in  space  relative  to  another  arbitrary  point.  While  this  may  sound 
simple,  it  is  extremely  important  because  it  allows  for  the  comparison  of  two  objects 
within  the  same  reference  frame.  This  comparison  and  the  forces  which  act  between  and 
on  the  objects  is  normally  the  true  goal  of  the  investigation. 

Typically  a  reference  frame  and  coordinate  system  are  chosen  in  such  a  way  as  to 
simplify  the  problem.  An  example  which  is  directly  applicable  to  this  thesis  is  the 
relative  reference  frame.  In  comparing  two  objects  which  are  both  in  motion  relative  to 
an  inertial  reference  frame  it  greatly  simplifies  the  problem  if  a  reference  frame  and 
coordinate  system  is  fixed  to  one  object. 

For  this  thesis  three  different  reference  frames/coordinate  systems  are  applicable. 
Primarily  among  these  is  the  inertial  reference  frame  in  which  all  physical  laws  are  based. 
Since  the  two  objects  of  interest  are  both  orbiting  the  Earth  it  is  natural  to  remove  the 
rotation  of  the  Earth  about  its  axis  as  a  simplification.  This  involves  an  Earth  Centered 
Earth  Fixed  (ECEF)  reference  frame/coordinate  system.  The  relative  reference  frame 
uses  one  of  the  spacecraft  as  the  center  of  the  frame  in  order  to  eliminate  any  identical 
orbital  motion  between  the  two  crafts.  In  all  three  cases  transformation  have  been 
developed  in  order  to  translate  between  any  of  the  coordinate  systems  without  loss  of 
accuracy. 
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A.  THE  M-50  INERTIAL  REFERENCE  FRAME 


Inertial  by  definition  is  a  reference  frame  which  is  fixed  in  space  or  is  moving 
with  a  constant  velocity.  The  Inertial  Reference  Frame  used  by  this  thesis  was  the  NASA 
standard  M-50  frame.  The  choice  of  this  frame  was  obvious  because  the  data  supplied  by 
NASA  for  this  thesis  was  given  in  terms  of  these  coordinates.  It  is  displayed  in  Figure 
2.1,  as  a  three  dimensional,  rectilinear  system  defined  by  x,  y  and  z  axes.  The  x  axis  is 
defined  by  the  Earth’s  Vernal  Equinox  or  the  First  Point  of  Aries.  The  z  axis  coincides 
with  the  Earth’s  rotation  vector,  and  roughly  points  toward  the  star  Polaris.  The  y  axis  is 
defined  by  the  cross  products  of  the  x  and  z  axis  unit  vectors. 


Polaris 


Figure  2.1  Inertial  Reference  Frame 

The  M-50  reference  frame  is  Earth  centered  but  the  Earth  does  rotate  relative  to 
the  frame.  Because  the  Earth  also  is  in  motion  about  the  Sun  it  could  be  argued  that  this 
is  not  truly  an  inertial  reference  frame.  While  this  is  true,  a  Geocentric  Inertial  Frame  is 
sufficiently  accurate  for  this  investigation  as  long  as  the  speeds  involved  are  not  very 
large  relative  to  the  speed  of  light,  because  then  that  Newtonian  mechanics  must  be 
replaced  with  Einsteinian. 

B.  WGS-84  EARTH  CENTERED  EARTH  FIXED  REFERENCE  FRAME 

In  order  to  eliminate  unnecessary  complications  caused  by  the  Earth’s  rotation  an 
Earth  Centered  Earth  Fixed  reference  frame  is  used.  The  use  of  this  reference  frame  was 
also  necessitated  by  the  fact  that  the  GEM-9  Gravity  Model  uses  the  WGS-84  coordinate 
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system.  Again,  as  seen  in  Figure  2.2,  the  WGS-84  is  a  three  dimensional,  rectilinear 
system  defined  by  x,  y  and  z  axes.  The  z  axis  is  defined  by  the  Earth’s  spin  axis  while  the 
x  axis  is  in  this  case  defined  by  0°  Latitude  and  Longitude.  Again  the  y  axis  is  defined  by 
the  cross  product  of  the  x  and  z  unit  vectors.  The  entire  reference  frame  is  in  rotation 
about  the  z  axis  at  a  rate  equal  to  the  Earth’s  rotation  rate. 


North  Pole 


Figure  2.2  Earth  Centered,  Earth  Fixed  Coordinate  System 

C.  RELATIVE  MOTION  REFERENCE  FRAME 

The  relative  motion  reference  frame  is  where  most  of  the  work  takes  place. 
Because  we  are  looking  at  two  objects  orbiting  the  Earth  in  near  identical  paths,  much  of 
the  orbital  motion  can  be  ignored.  It  is  logical  to  assume  that  forces  acting  on  the  two 
objects  identically  in  an  inertial  or  ECEF  reference  frame  can  also  be  ignored  with  little 
loss  in  accuracy  when  a  relative  system  is  adopted.  Whether  this  assumption  is  true  or 
not  is  the  subject  of  this  thesis. 

There  are  two  potential  coordinate  systems.  The  first  is  called  Local  Vertical, 
Local  Horizontal  (LVLH)  and  is  defined  relative  to  the  orbit  of  the  spacecraft,  while  the 
second  is  called  Roll,  Pitch  and  Yaw  (RPY),  is  defined  relative  to  the  spacecraft  itself. 
Both  have  advantages  or  disadvantages  depending  on  what  is  required. 

1.  Local  Vertical.  Local  Horizontal  Coordinate  System 

LVLH  is  again  a  rectilinear  system  with  x,  y,  and  z  coordinates.  The  x  axis  is 
defined  by  the  direction  of  motion  while  the  z  axis  points  radially  outward  from  the  Earth. 
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The  y  axis  is  defined  by  the  cross  product  of  the  x  and  z  axis.  As  seen  in  Figure  2.3  the 
horizontal  plane  is  defined  by  the  x  and  y  axis  and  the  vertical  plane  by  the  x  and  z  axis. 
A  LVLH  coordinate  system  is  in  constant  rotation  relative  to  an  ECEF  coordinate  system 
(typically  about  the  y  axis.)  The  magnitude  of  that  rotation  is  equal  to  the  spacecraft’s 
angular  velocity  about  the  Earth’s  center.  Since  the  LVLH  coordinate  system  is  based 
upon  the  object’s  orbit  and  not  it’s  attitude  and  because  the  orbiting  objects  are  small 
relative  to  the  earth,  they  can  normally  be  treated  as  point  masses. 

y 


Figure  2.3  Local  Vertical,  Local  Horizontal  Coordinate  System 

2.  Roll,  Pitch  and  Yaw  Coordinate  System 

The  RPY  coordinate  system  is  identical  to  that  of  an  aircraft’s  with  the  x  or  roll 
axis  pointing  out  the  nose  of  the  spacecraft,  the  z  or  yaw  axis  out  the  spacecraft’s  belly 
and  the  y  or  pitch  axis  pointing  out  the  starboard  wing.  See  Figure  2.4  for  a  graphical 
representation  of  this  coordinate  system.  This  is  a  very  advantageous  system  for  any 
sensors  or  personnel  located  on  the  spacecraft,  but  since  the  spacecraft  may  be  in  rotation 
about  its  center  of  mass  it  is  not  ideal  when  discussing  the  motion  between  two  orbiting 
objects. 
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Figure  2.4  Roll,  Pitch  and  Yaw  Coordinate  System 


D.  COORDINATE  TRANSFORMATIONS 

As  discussed  previously,  it  is  critically  important  to  be  able  to  translate 
coordinates  and  forces  between  the  various  reference  frames  in  order  it  insure  the 
accuracy  of  the  Equations  of  Motion  used.  The  translation  process  can  be  a  complex  one 
involving  detailed  math,  resulting  in  the  addition  of  many  new  forcing  terms  into  the 
Equations  of  Motion.  These  pseudo-forces,  while  having  no  existence  in  inertial  space 
are  very  important  in  a  rotating  reference  frame.  Fortunately,  many  of  these  terms  cancel 
when  dealing  with  relative  motion. 

The  fact  that  all  coordinate  systems  used  are  rectilinear  is  convenient  because 
transformation  matrices  between  systems  exist  and  are  well  known.  Typically  a 
transformation  can  be  view  as  a  series  of  three  rotations  about  the  coordinate  axis  such 
that  one  coordinate  system  is  translated  into  the  other.  The  rotation  about  a  single  axis 
can  be  seen  in  Figure  2.5  . 
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a. 


Figure  2.5  Axis  Rotation 


Mathematically  this  rotation  can  be  done  by  multiplying  a  vector  by  a  rotational 
matrix  to  get  the  rotated  coordinates.  This  rotational  matrix  is  determined  from  the  unit 
vectors  of  the  two  axis. 

al-bl  ax  •  b2  ax-b2  cos  <j)  sin0  0 

ACB=  a2-bx  a2-b2  a2-b3  =  -sin0  cos  (j)  0  (2.1) 

a2  -bx  a3-  b2  a3  ■  b2  0  0  1 

The  rotational  matrix  in  Equation  2.1  can  then  be  multiplied  with  two  other 
rotational  matrices  to  get  the  final  transformation.  The  order  of  the  multiplication  is 
important  as  are  the  order  of  the  rotations,  to  insure  that  the  x  axis  in  the  first  coordinate 
system  lines  up  with  the  x  axis  in  the  new  coordinate  system. 

‘CE  =[,CA\ACB\BCE]  (2.2) 

Transformation  matrices  were  generated  and  used  in  translations  between  all 
coordinate  systems  to  insure  proper  accuracy  of  the  data. 
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III.  THE  PROPAGATORS 


The  purpose  of  any  propagator  is  to  accurately  predict  a  spacecraft’s  position  and  velocity 
at  some  future  point  in  time.  There  are  many  different  mathematical  methods  to  produce  these 
future  vectors,  each  with  their  own  advantages  and  shortcomings.  While  they  may  differ  in 
application  or  procedure,  all  claim  Newton's  Universal  Law  of  Gravity,  Equation  3.1,  as  their 
source. 


r,  Gmim2 
F  = - — 

r2  (3.1) 

where  G  is  the  Universal  Gravitational  Constant  and  r  is  the  distance  separating  the  point  masses 
mi  and  m2. 

To  properly  apply  this  physical  law  to  the  case  of  orbital  mechanics,  some  general 
assumptions  must  be  made.  The  first  is  that  the  two  objects  of  interest  are  point  masses  which  is 
valid  if  we  assume  the  bodies  are  either  spherically  symmetrical  or  relatively  small.  Secondly,  it 
is  assumed  that  these  two  point  masses  are  the  only  ones  in  the  system.  In  the  case  of  an  object 
orbiting  the  Earth,  the  assumption  of  two  point  masses  gives  a  sufficiently  accurate  answer. 
Given  these  underlying  assumptions  and  applying  Newton's  Third  Law  of  Motion  we  obtain 
Equation  3.2  which  is  the  equation  of  the  force  on  an  orbiting  spacecraft. 


mr  =  — 


GMm  _ 


(3.2) 


Dividing  out  the  spacecraft  mass  we  have  Equation  3.3  which  governs  motion  for  the 
Restricted  Two  Body  Problem.  The  restriction  is  due  to  the  presupposed  assumption  that  the 
Earth’s  mass  M  »  m  which  in  the  case  of  a  spacecraft  orbiting  the  Earth  is  definitely  true.  This 
assumption  additionally  fixes  the  coordinate  system  at  the  Earth’s  center  since  the  spacecraft  will 
have  very  little  effect  of  the  system’ s  Center  of  Mass. 
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(3.3) 


n  GM  _ 
r  +  — —r  =  0 


The  zero  to  the  right  of  the  equal  sign  in  Equation  3.3  results  from  our  assumptions  that 
there  are  only  two  bodies  in  the  system  providing  accelerating  forces.  This  equation, 
unfortunately,  only  approximately  describes  the  motion  of  a  spacecraft  about  the  Earth.  In  reality 
there  are  small  perturbations  caused  by  other  planets,  atmospheric  drag,  the  fact  that  the  Earth  is 
not  really  a  point  mass,  etc.  An  expression  of  the  form, 


r  =  ■ 


GM  _ 

— rr+ap 


(3.4) 


where  ap  is  the  sum  of  all  the  perturbing  accelerations,  is  more  accurate.  Whether  to  include 
certain  accelerations  in  ap  depends  on  the  order  of  magnitude  of  these  perturbing  forces  acting  on 
the  body,  and  on  the  time  interval  over  which  they  are  allowed  to  affect  its  motion. 

There  are  a  plethora  of  possible  solution  methods  used  to  solve  the  second  order 
differential  Equation  3.4  (Ref.  2).  The  distinguishing  features  between  the  different  propagators 
is  the  manner  in  which  this  equation  is  solved.  For  the  purposed  of  this  thesis  only  three 
propagators  are  compared,  they  being  the  most  popular. 

A.  THE  COWELL  PROPAGATOR 

The  Cowell  propagator  is  the  highest  fidelity  propagator  used.  It’s  accuracy  is  due  to  the 
fact  that  it  solves  the  complete  unsimplifed  differential  equation.  It  begins  by  breaking  the 
second  order  differentials  into  two  coupled  first  order  different  equations. 


r  =  v 


v  =  — 


GM 


,-r+ar 


(3.5) 


These  equations  are  solved  simultaneously  using  a  Runge-Kutta  4th  order  integrator.  A  built  in 
feature  of  the  code  decreases  the  time  step  used  in  the  numerical  integrator  scheme  when  the 
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solution  varies  to  rapidly  over  a  short  time  interval.  The  step  size  is  chosen  in  an  effort  to 
contain  numerical  truncation  and  rounding  errors  below  a  specified  tolerance  level. 

While  the  Cowell  propagator  is  generally  the  most  accurate  it  is  also  the  slowest.  This  is 
because  it  must  perform  the  numerical  computations  over  a  relatively  large  number  of  small  time 
steps  in  order  to  maintain  accuracy.  In  addition,  the  time  required  significantly  increases  when 
atmospheric  and/or  gravity  models  are  used  to  increase  the  propagator’s  accuracy. 

B.  THE  CLOHESSY- WILTSHIRE  PROPAGATOR 

This  propagator  uses  a  set  of  linearized  Equations  of  Motion  for  relative  motion 
developed  by  Clohessy-Wiltshire.  The  derivation  is  taken  from  Reference  2  and  Reference  3. 
The  solution  of  the  equations  provides  the  position  and  velocity  of  a  spacecraft,  the  chaser,  with 
respects  to  another,  the  target.  In  order  to  define  vectors  between  these  two  objects  an 
orthogonol  reference  frame  is  attached  to  the  target  and  moves  with  it.  The  coordinate  system  is 
a  form  of  Local  Vertical,  Local  Horizontal. 

The  vector  position  of  the  chaser  is  given  by 

r  =  rT+p  (3.6) 

which  is  differentiated  with  respects  to  the  M50  inertial  coordinate  system  giving 

?  =  fr+p+2(  coxp)+5)Xp+  ebx(  obxp)  (3.7) 

where  (omega)  is  the  angular  velocity  of  the  target  about  the  center  of  the  Earth.  If  we  take 

r  =  g  +  ap  (3.8) 

where  g  is  the  acceleration  due  to  gravity  and  A  is  the  acceleration  due  to  external  forces, 
Equations  3.7  and  3.8  can  be  resolved  into  x,  y  and  z  components. 
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(3.9) 


x  =  —g—  +  Ax  +  2a>z  +  o)z  +  g)2x 
r 

y  =  -g-+Ay 
r 

( z  +  .  , 

-g -  +Az  +  gT-2cax-CQx  +  cQz 

V  r 


Assuming  that  the  distance  between  the  target  and  chaser  is  much  smaller  than  the  orbital  radius 
of  the  target,  i.e  (p)2 « r2j ,  the  following  approximations  can  be  made 


r-rT[1+-rj 

g-gjl-— 

V  rrj 

x  x 

~g~  ®  ~gT— 
r  rT 

y  y 

g  ~  ~gT 
r  rr 


Z  +  rr 


Sr  1  — 


(3.10) 


which  by  substitution  into  Equation  3.9  yields 


x  =  -g — I -  ax+2coz  +  coz  +  co2  x 
rr 


y  gr  Ay 
rr 


(3.11) 


z=2  gT — Az  +  2(ax  +  <m  +  co2z 
rr 
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If  we  make  the  additional  assumption  that  the  target's  orbit  is  circular,  then 

(0=0 
(0= 

and  Equation  3.11  becomes 

x  =  Ax  +  2caz 

y  =  Ay-(02y  (3-12) 

z  =  Az-2cok  +  3co2z 

Since  y  is  uncoupled  from  x  and  z  in  Equation  3.12,  it  can  be  solved  simply  and  separately. 
Assuming  a  solution  of  the  form 

y  =  Asintot  +  Bcosox+^- 

of  (3.13) 

A  and  B  are  determined  from  the  initial  conditions,  giving  the  solution 

y 

y  =  y  cosCOt  +  —  (Osin(Ot-\ — ffi-costuf ) 

CO  CO2  (3.13) 

A  • 

y  =  yocos(Ot  -  yj  osintut )+  — ^-sinfflt 

CO 

The  equations  for  x  and  z  are  more  difficult  because  they  are  coupled.  We  begin 
integrating  the  equation  for  x  once  and  substitute  the  result  into  the  equation  for  z. 

z  +  2(o(2coz  +  AJ  +  cj)-3(o2z-Az  =  0  (3.15) 

Assuming  that  the  solution  to  Equation  3. 15  is  of  the  form 

z-  Asmcot  +  Bcoscot - - - —  H — -  (3.16) 

(O  (0(0 
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it  can  be  substituted  into  the  equation  for  x  dot  resulting  in: 


x  =  2  coA  sin  cot  +  2  coB  cos  cot  -  3  Ax  t  -  3  a  +  2  - 


(3.17) 


Integrating, 


x  —  -2A  cos  ax  +  lBsinax  —  Axt2  —  3  at  +  2 — t  +  ci 

2  co 


(3.  IB) 


and  applying  the  initial  conditions  we  obtain 


1  2 

A  =  —(  Zo+—Ax) 

CO  CO 

r,  _2_  .  2 _ A± 

CO  of 

'  Ci  *■“  Xo  ^  co 

2  2 

Ci  -  x<A — ( Zo H —  Ax) 
CO  CO 


(3.19) 


Substituting  and  reordering  terms, 


*=*  +z  [6(or— sinor)l+4  3*  +z9  -(1-cosor)  +4  ~r  +4r(l-cosor)  +4  — (or-sinor) 

L  \_  co  j  \co  J  'L  2  co~  J  ico 

r  -i  .  fsincurl  A,  , 

y=y  [cosor]+;y0 -  +— (1 -coscar) 

L  co  J  cor 

z  =  k  -—(1-cos cct)  +z J4-3cosar]+4  +4  — {ca -sin cot)  +A  -(1-cosor) 

°]_  co  J  a  J  l  car  J  Lor 

r 4  i  \2  i 

i: = 4  [4  cos  or  -  3]  -f-  4  [6o)(  1  -  cos  or)] + 4[2  sin  or]  +  4  sin  or  -  3f  +4  -(1-cosor)  (3.20) 


L  \ 

v=y0cosor-y0(o)sinor)  +— sin  or 

co 


r  f  2,  ."I  fsinor 

z  =  x  i-2sinorl  +  z  f3o)sinor]-r4Lcosor]  +  4| — (1-cosor)  -f  A - 

"L  J  oL  -  °  L  co  J  L  or  J 
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The  Clohessy-Wiltshire  equations  are  convenient  since  they  give  answers  in  LVLH 
coordinates.  This  implies,  assuming  the  radius  and  velocity  vectors  are  obtained  by  the  chasing 
vehicle,  that  no  coordinate  transformation  is  required  for  orbit  determination.  A  drawback  the  to 
the  solution  given  above  is  that  they  are  based  upon  the  assumption  that  the  target  orbit  is 
circular.  This  produces  an  inherent  error  in  the  algorithm  whose  magnitude  varies  with  distance 
between  the  vehicles  and  eccentricity  of  the  two  orbits.  The  greater  the  distance  and  eccentricity 
the  greater  the  magnitude.  These  disadvantages  most  likely  make  the  Clohessy-Wiltshire 
Propagator  the  least  accurate  of  the  propagators  used. 

C.  THE  ENCKE  PROPAGATOR 

While  the  Cowell  method  is  probably  the  most  accurate  propagation  method  it  can  be 
cumbersome  and  time  consuming,  especially  if  the  perturbing  acceleration  is  much  smaller  than 
the  force  of  gravity.  Because  this  is  the  case  for  an  Earth  orbiting  spacecraft  an  alternative 
technique  based  upon  integrating  the  vector  difference  between  the  Two  Body  Solution  and  the 
perturbed  solution  rather  than  the  explicit  state  vector  found  in  the  Cowell  Method.  This  is  the 
idea  behind  the  Encke  Method.  The  particular  version  of  the  Encke  method  developed  here  is 
taken  from  Reference  4. 

The  Encke  Method  assumes  that  the  spacecraft  is  traveling  in  a  conical  path  called 
osculating  orbit.  This  osculating  orbit  is  determined  by  the  position  and  velocity  vectors  at  a 
given  time  fusing  the  classical  Two  Body  Problem.  Because  of  the  perturbing  accelerations,  the 
actual  orbit  diverges  in  time  from  the  osculating  orbit  and  these  deviations  are  what  is  found  by 
the  method.  These  vector  differences  are  then  added  to  the  osculating  orbit  position  and  velocity 
vectors  to  give  the  actual  state  vectors  at  a  future  time. 

As  stated  previously,  the  osculating  orbit  is  derived  using  the  restricted  two  body.  The 
method  employed  in  Reference  4  is  based  upon  the  so-called  “f  ’  and  “g”  functions.  These  are 
determined  from  the  position  and  velocity  vectors  at  a  prescribed  time  t0,  and  can  be  used  to 
generate  the  classical  orbital  element.  The  f  and  g  functions  are  given  in  Equation  3  (Ref.  4). 
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/  =  -[(cos(E-£„)-l]  +  l 


=  J^[sin(E-E,)-(E-Ej]  +  (t-0 


f  =  -^sin(E-EJ 

ror 

8  =  !--[!-  Cos(£-£0)] 


(3.21) 


Position  and  velocity  vectors  of  the  osculating  orbit  are  given  by: 


r„c  =  fre  +  gv0 

Ksc  =  fro+fr0 


(3.22) 


In  the  Encke  formalism,  these  vectors  are  used  in  the  differential  equation  for  the 
difference  vector  (8): 


d2S  /it  /I  (,  rM_ 
dt2  r3  r3  {  r3  I 

***  ' osc  osc  \  / 


(3.23) 


To  avoid  difficulties  in  evaluating  r,  slight  modifications  to  Equation  3.23  are  made. 


Since 


r(t)  =  rosc(t)  +  S(t) 


(3.24) 


we  have 


1— f-  =  -F(*)=l-(l  +  *)’ 


(3.25) 


where 


8  •  -  27 ) 


(3.26) 
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The  function  F(q)  can  then  substituted  into  Equation  3.23.  This  revised  second  order  differential 
equation  can  now  be  integrated  to  find  8  and  its  first  derivative  v.  When  these  are  substituted 
into  Equation  3.27,  the  position  and  velocity  vectors  for  the  desired  time  are  found: 


r  =  r  +  8 


v  =  vOJC  +  v 


(3.27) 


While  this  propagator  is  more  accurate  than  the  Clohessy-Wiltshire  Propagator  since  it 
does  not  depend  upon  to  circular  orbits  it  does  restrict  itself,  due  to  the  method  of  determining 
the  osculating  orbit,  to  planar  orbits.  This  means  that  greater  error  is  likely  in  the  out  of  plane 
direction,  the  y  axis  in  the  LVLH  reference  frame.  Whether  this  error  is  significant  is  yet  to  be 
determined. 
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IV.  THE  PERTURBING  ACCELERATIONS 


The  perturbing  accelerations  which  affect  the  motion  of  a  spacecraft  in  Earth’s 
orbit  are  many.  They  are  caused  by  physical  phenomenon  existing  outside  the  orbital 
system,  such  as  solar  radiation,  because  as  well  as  due  to  the  fact  that  many  of  the 
assumptions  made  in  developing  the  Restricted  Two  Body  Problem  are  too  simplistic. 
The  simplifying  assumptions  inherent  in  the  Two  Body  Problem,  such  as  disregarding  the 
existence  of  three  body  effects,  cause  errors  in  the  predictions.  Regardless  of  the  source, 
these  perturbing  accelerations  should  be  considered  in  the  propagation  techniques  in  order 
to  insure  accurate  predictions. 

A.  THE  PERTURBATIONS 

When  considering  a  Low  Earth  Orbit  (LEO),  just  about  every  conceivable 
perturbing  force  has  an  effect.  Each  of  them  must  be  reviewed  and  either  incorporated 
into  the  propagator  or  discarded  depending  upon  the  magnitude  of  the  perturbation  they 
produce. 

1.  Solar  Pressure. 

The  perturbation  known  as  solar  pressure  is  caused  by  the  impingement  of  solar 
photons  on  the  surface  of  a  spacecraft.  One’s  first  impulse  is  to  discard  this  as  negligible. 
Before  doing  so,  one  should  consider  that  in  space  there  is  no  counter  force,  such  as 
atmospheric  drag,  to  resist  such  a  force.  Over  a  sufficiently  long  time  interval  this 
perturbing  force  could  build  up  to  a  significant  level.  Near  Earth,  the  Solar  Pressure 
Constant  is  4.7  x  10'5  dynes/cm2. 

Typically  Solar  pressure  effects  the  satellite’s  orbital  eccentricity,  but  as  implied 
above,  it  typically  takes  time  for  it’s  effects  to  be  felt.  The  time  frame  over  which  the 
propagator  must  be  accurate  for  the  purpose  of  this  research  is  one  orbit,  or  approximately 
90  minutes,  and  therefore  this  effect  can  be  reasonably  neglected. 
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2.  Magnetic 


The  Earth’s  magnetic  field,  like  any  magnetic  field,  can  induce  a  force  on  a 
magnetized  or  a  current  carrying  body.  The  magnitude  of  this  force  is  dependent  upon 
the  size  of  the  magnetic  fields  and/or  current.  Such  a  force  is  often  used  by  spacecraft  via 
magnetic  torquers  to  actively  introduce  a  change  in  the  spacecraft  attitude.  Even  if  there 
were  no  active  use  of  the  magnetic  field,  a  spacecraft,  after  a  sufficient  time,  develops  a 
slight  magnetic  field. 

The  difficulties  in  using  the  Earth’s  magnetic  field  for  atitude  control  is  really  a 
boon  for  an  orbital  propagator  because  the  constantly  changing  direction  and  magnitude 
which  creates  great  difficulties  in  developing  a  lengthy,  direction-stable  force  does  not 
create  a  very  coherent  perturbing  force.  Additionally,  because  the  Earth’s  magnetic  field 
is  so  weak,  only  a  fairly  large  current  in  the  orbiter  creates  an  appreciable  effect.  For 
these  reasons,  magnetic  perturbations  are  ignored  in  this  thesis. 

3.  Third  Body  Effects 

The  reality  of  orbital  mechanics  is  that  there  are  several  gravitating  bodies  rather 
that  the  two  assumed  in  the  Restricted  Two  Body  Problem.  The  Moon,  the  Sun  and  the 
other  planets  all  have  a  gravitational  effect  on  an  orbiting  spacecraft.  To  accurately 
model  these  effects  the  system  should  be  expanded  from  two  bodies  to  three,  four,  etc. 
This  was  attempted  for  a  three  body  system  and  was  found  to  be  analytically  intractable 
except  in  very  limited  circumstances.  Extensions  to  a  system  of  four  or  greater  bodies 
turns  out  to  be  of  very  limited  value  for  practical  matters  involving  orbiting  satellites. 

One  possibility  is  to  treat  these  effects  like  any  other  perturbation.  It  has  been 
found  that  the  accumulated  effects,  while  they  do  exist,  are  again  insignificant  over  the 
short  duration  of  a  single  orbit.  Therefore  these  perturbations  are  also  ignored. 


22 


4.  Earth  Gravity  Harmonics 


The  Restricted  Two  Body  Problem  assumes  that  the  Earth  is  spherically 
symmetrical  and  therefore  can  be  represented  by  a  point  mass.  In  reality  it  is  a  large 
semi-homogenous  ball.  For  a  point  mass,  it  can  be  assumed  that  the  force  of  gravity  is 
constant  at  every  point  on  a  sphere  encompassing  that  point  mass.  For  a  semi- 
homogenous  ball,  every  point  on  that  same  sphere  can  experience  a  different  gravitational 
force.  In  order  to  take  this  variation  into  account  a  model  of  this  variation  has  been 
developed.  In  the  case  of  gravity,  the  concept  of  a  gravitational  potential  provides  this 
model. 


By  mapping  the  gravitational  potential  at  all  points  around  the  Earth  it  is  possible 
to  determine  the  acceleration  due  to  gravity  at  any  such  point.  The  acceleration  in  a  given 
direction  is  simply  the  gradient  of  the  potential  in  that  particular  direction.  The  external 
geopotential  is  based  upon 


where: 


1  +  X  S  f  ■ -]  Pnm  (Sin  4) (Cnm  C°S  +  $nm  «*) 
n= 2  m= O'-  ^  ' 


(4.1) 


|i  =  The  Earth’s  Gravitational  constant 

r  =  The  magnitude  of  the  spacecraft  radius  vector 

a  =  The  Semi-major  axis  of  the  Central  ellipsoid  (Earth) 

n,m  =  Degree  and  order  of  each  term 

cj)  =  Geocentric  latitude 

%  =  Geocentric  longitude 

Cnm,  Snm  =  Gravitational  coefficients 

Pnm  =  Associated  Legendre  functions 

N  =  number  of  terms  used 
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The  Gravitational  Coefficients  can  be  broken  into  two  categories,  Cn0  ,  the  zonal 
harmonics  or  the  Cnm  ,  the  tesseral  and  sectorial  harmonics.  The  zonal  harmonics  are 
also  know  as  the  Jn  coefficients  of  which  J2  is  the  most  dominant.  The  harmonics  define 
the  shape  and  distribution  of  the  Earth  mass. 

The  J2  effect  can  have  a  significant  effect  on  a  spacecraft’s  orbit  in  a  short  amount 
of  time.  For  this  reason  a  gravity  model  must  be  examined  to  see  if  it  has  any  substantial 
effect.  Also  it  must  be  determined  to  what  fidelity,  i.e.  number  of  coefficients,  the  model 
should  have. 

5.  Atmospheric  Drag 

The  final  perturbing  force  considered,  is  that  due  to  atmospheric  drag.  Even  as 
high  as  1000  km  the  Earth  has  a  measurable  atmosphere.  While  it  is  true  that  it  is 
essentially  only  an  indistinct  plasma,  there  are  still  particles  which  can  impact  a 
spacecraft  and  cause  a  change  in  velocity.  Also,  as  the  altitude  decreases  the  atmosphere 
becomes  more  dense  and  therefore  more  of  a  factor. 

For  LEO  spacecraft,  atmospheric  drag  can  be  the  dominant  perturbing  force.  The 
drag  a  spacecraft  feels  is  a  function  of  velocity,  atmospheric  density  and  the  spacecraft’s 
ballistic  coefficient.  The  velocity  is  constant  for  a  circular  orbit,  but  as  discussed 
previously,  orbits  are  never  truly  circular.  For  elliptical  orbits  velocity  is  constantly 
changing  depending  on  the  point  along  the  ellipse.  The  fact  that  the  atmosphere’s  density 
is  also  variable  as  a  function  of  both  position  and  time  further  complicates  the  problem. 
The  ballistic  coefficient  can  be  assumed  constant  since  it  depends  upon  the  shape  of  the 
Orbiter  which  does  not  change. 

The  fact  that  velocity  and  density  are  variable  means  that  drag  is  constantly 
changing,  typically  in  a  very  short  time.  As  a  result  of  this  rapid  change,  Atmospheric 
Drag  has  the  potential  to  significantly  effect  an  orbital  prediction.  Velocity  is  easily 
gained  from  an  osculating  orbit  model,  but  density  often  requires  a  very  complex  model. 
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B.  THE  MODELS 


For  the  purpose  of  this  thesis  only  two  perturbing  forces  are  used,  Earth 
Gravitational  Harmonics  and  Atmospheric  Drag.  There  are  numerous  models,  each  with 
their  own  advantages  and  disadvantages.  Models  for  this  thesis  were  primarily  chosen 
because  of  their  availability,  but  fortunately,  both  are  also  extremely  accurate. 

1.  The  GEM-9  Geopotential  Model 

All  Geopotential  Models  are  essentially  the  same.  Equation  4.1  is  used  to  generate 
the  gravitational  potential  at  a  specific  location.  They  differ  in  the  values  of  the 
gravitational  coefficients  and  the  number  used.  For  this  thesis  the  GEM-9  model  using 
coefficients  from  Table  3  of  Reference  5  are  used.  The  model  has  a  maximum  resolution 
of  30  x  30,  and  uses  the  WGS-84  coordinate  system.  This  is  an  ECEF  standard 
coordinate  system. 

2.  The  Jacchia  71  Atmospheric  Model 

There  are  a  number  of  potential  Atmospheric  models.  Some  are  purely 
exponential,  while  others  are  more  complex,  taking  into  account  the  date,  time,  latitude 
and  longitude  and  the  effect  of  the  sun  on  the  upper  atmosphere.  The  Jacchia  Models  are 
of  the  latter  category  making  them  large  but  more  accurate.  This  model  was  developed  by 
Jacchia  in  1970  and  takes  into  account  the  11 -year  and  27-year  solar  cycles.  It  has  the 
added  advantage  of  using  the  WGS-84  coordinate  system. 
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Y.  COMPARING  PROPAGATORS 


A.  THE  DATA  SET 

The  easiest  way  to  judge  the  accuracy  of  any  algorithm  is  to  compare  it  to  actual 
data  from  the  phenomenon  being  modeled.  If  the  algorithm  results  in  a  relatively  good 
match  with  the  data  the  model  and  algorithm  are  considered  acceptable.  The  usual  tack 
taken  by  a  modeler,  is  to  start  with  the  simplest  model  available  to  describe  a 
phenomenon,  and  include  effects  (or  fine  tune  the  model)  until  such  a  time  as  the  model 
is  acceptable.  It  may  happen  that  further  refinement  is  possible  but  perhaps  unnecessary 
to  produce  a  workable  model  and  algorithm.  If  not,  it  is  generally  back  to  the  drawing 
board  for  the  individual  who  is  developing  the  model  and  algorithm. 

The  data  set  in  this  thesis  comes  from  a  rendezvous  mission  performed  by  the 
Space  Shuttle,  and  consists  of  a  series  of  state  vectors  for  the  Space  Shuttle  and  a 
rendezvous  target.  Each  state  vector  is  comprised  of  the  spacecraft’s  orbital  position  and 
velocity  as  well  as  a  time  stamp  when  the  measurement  was  made.  The  entire  data  set  has 
3238  state  vectors  taken  at  approximately  3.5  seconds  intervals.  The  data  set,  as  seen  in 
Figure  5.1,  contains  a  record  of  two  orbits  during  which  time  the  Orbiter  performs  a 
rendezvous  with  the  target  vehicle.  Unfortunately,  it  was  determined  near  the  end  of  the 
project  that  the  data  could  have  been  processed  using  a  high-fidelity  propagator  to  fill  in 
the  gaps.  This  fact  limits  the  possible  accuracy  to  that  provide  by  the  high-fidelity 
propagator.  It  also  could  produce  a  more  accurate  prediction  to  appear  less  accurate. 

In  Figure  5.1,  the  data  is  presented  as  x  versus  z,  but  since  the  x  coordinate 
closure  rate  is  constant,  the  figure  also  demonstrates  orbital  positions  as  a  function  of 
time.  All  state  vectors  are  provided  in  M50  inertial  coordinates.  Additionally,  because 
the  data  is  of  an  actual  rendezvous,  there  are  numerous  thruster  firings  performed  by  the 
Orbiter  during  the  time  interval  covered  by  the  data  set.  These  firings  produce  anticipated 
errors  in  the  propagators,  but  since  they  are  essentially  impulses  they  are  also  readily 
identifiable. 
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Figure  5.1  Data  Set  in  LVLH  Coordinates 


B.  THE  PROGRAMS 

Since  the  motivation  for  this  thesis  is  the  efficacy  of  various  propagators  in  the 
development  of  the  RPOP  software,  it  seemed  logical  to  retain  the  C++  coding  of  the 
actual  candidates.  The  majority  of  the  code  used  in  this  thesis  is  the  NPSLite  software 
developed  by  LT  Lee  Barker  (Ref.  6),  but  numerous  modifications  and  additional  code 
has  been  written  by  the  author. 

The  Cowell  propagator  used  was  originally  developed  by  Robert  Gottlieb  and 
Mike  Fraietta  of  McDonnell  Douglas  Corporation  for  NASA.  Some  code  modification 
were  made  by  LT  Lee  Barker  (Ref.  6)  to  make  the  algorithm  work.  Additional 
modifications  were  incorporated  by  the  author  into  the  C++  code  written  for  this  thesis, 
but  no  changes  to  the  core  algorithm  were  made.  In  order  to  properly  compare  the 
performance  of  the  Cowell  propagator  with  the  Clohessy-Wiltshire  propagator,  both  the 
orbiter  and  the  target  state  vectors  were  predicted  and  converted  to  the  LVLH  reference 
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frame.  This  was  to  nullify  any  common  errors  made  by  the  algorithms  which  might  have 
occurred  in  both  spacecraft  predictions 

The  Clohessy-Wiltshire  propagator  was  originally  written  by  Scott  Tamblyn  of 
NASA’s  Johnson  Space  Center  and  Jack  Brazzel  of  McDonnell  Douglas.  Modifications 
to  this  algorithm  were  made  to  incorporate  perturbing  accelerations.  The  perturbing 
accelerations  include  a  Jacchia  70  atmospheric  model  and  a  GEM-9  gravity  model  coded 
for  use  by  the  Cowell  propagator.  By  using  the  same  perturbing  accelerations  as  the 
Cowell,  a  potential  source  of  error  was  eliminated.  Both  the  orbiter  state  vector  and  the 
perturbing  accelerations  had  to  be  converted  from  the  inertial  frame  to  the  LVLH  frame 
prior  to  use  by  the  Clohessy-Wiltshire  propagator. 

The  Encke  Propagator  used  the  “f  and  g”  functions  developed  by  LT  Les 
Makepeace  (Ref.  1)  and  coded  by  LT  Barker  (Ref.  6)  for  NPSLite.  A  Runge-Kutta  4-5 
integrator  which  integrates  Equation  3.23  was  incorporated  by  the  author.  Gravitational 
and  atmospheric  models  were  once  again  taken  from  the  Cowell  propagator  to  maintain 
coherence  in  the  project.  As  with  the  Cowell  propagator,  state  vectors  were  predicted  in 
inertial  coordinates  and  then  converted  to  LVLH  for  the  purpose  of  comparison. 

Numerous  difficulties  were  experienced  in  compiling  and  running  the  Encke 
propagator  when  the  gravitational  and  atmospheric  models  were  added.  Without  the 
additional  models  the  Encke  propagator  ran  well.  Addition  of  the  perturbing  forces 
caused  a  divergence  of  the  orbit  predictions.  Unfortunately  there  was  insufficient  time  to 
track  this  error  down,  and  therefore  an  interesting  portion  of  the  project  could  not  be 
accomplished. 

C.  PERFORMANCE 

The  performance  of  the  different  propagators  has  not  been  judged  solely  on  the 
accuracy  of  the  predicted  state  vectors,  but  also  on  how  that  accuracy  varied  as  the  time 
step  is  increased.  Additionally,  the  amount  of  computation  time  required  to  perform  the 
propagation  is  compared.  It  was  anticipated  that  in  order  to  obtain  the  most  accurate 
prediction,  a  high  fidelity  propagator  would  be  needed  to  break  up  the  prediction  time 
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interval  into  a  series  of  smaller  time  steps  over  which  each  propagator  could  be  run. 
Ideally,  the  propagator  requiring  the  shortest  computation  time  over  the  largest  acceptable 
time  step  would  optimize  performance. 

For  the  purpose  of  this  project  four  time  steps  were  chosen.  The  smallest  step  was 
based  on  the  separation  between  consecutive  state  vectors  in  the  data  set,  approximately 
3.84  seconds.  In  order  to  insure  that  actual  state  vectors  would  be  available  for 
comparison,  the  remaining  three  time  steps  were  simply  a  factor  two,  five  and  ten  times 
that  base  time  step.  This  resulted  in  time  steps  of  7.68,  19.2  and  38.4  seconds.  Average 
computation  time  was  determined  by  dividing  the  total  amount  of  time  required  to 
propagate  all  3238  state  vectors  over  all  four  time  steps  by  the  total  number  of  state 
vectors  propagated. 

D.  FINDINGS 

The  Cowell  and  Encke  propagators  proved  to  be  the  most  accurate  as  seen  by 
Figure  5.2  (a)  and  (b).  The  Encke  Coordinate  and  Velocity  Differentials  were  virtually 
identical  to  those  of  the  Cowell  method.  Differentials  is  defined  as  the  difference 
between  the  actual  state  vector  and  the  predicted  state  vector. 

The  times  represented  on  the  x-axis  of  the  two  figures  is  in  term  of  seconds  since 
the  beginning  of  the  calendar  year  accordance  with  the  time  step  of  each  state  vector. 


The  Ciohessy-Wiitshire  propagator  predicted  sinusoidal  variations  as  seen  in 
Figures  5.3  (a)  and  (b).  The  cause  of  the  sinusoidal  motion  is  a  result  of  the  initial 
assumption  in  deriving  the  Ciohessy-Wiitshire  Equations  that  the  target  orbit  is  circular. 
In  reality  the  target  orbit  is  more  nearly  elliptical,  thereby  introducing  an  error.  This  fact 
was  discussed  in  Chapter  III,  Section  B. 

Additionally,  large  errors  caused  by  thruster  firings  can  be  seen  in  both  Figures 
5.3  and  5.4.  The  largest  of  these  occur  at  identical  times  with  the  Encke,  Cowell  and 
Ciohessy-Wiitshire  propagators,  but  notice  that  the  Ciohessy-Wiitshire  is  more 
susceptible  to  small  firings  judging  form  the  relatively  large  number  of  spikes  occurring 
in  Figure  5.3  compared  to  those  in  Figure  5.2. 


Figure  5.3  (a)  Position  Differential  of  CW  Propagator  for  a  3.84  Time  Step 
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Figure  5.3  (b)  Velocity  Differential  of  CW  Propagator  for  3.84  Time  Step 
It  is  easier  to  view  the  separate  differences  in  the  three  propagators  when  they  are 
compared  with  respect  to  Position  and  Velocity  Components.  This  type  of  comparison 
also  makes  it  possible  to  see  the  differences  produced  by  adding  perturbing  accelerations 
to  the  Clohessy-Wiltshire  propagator.  For  the  purpose  of  this  analysis  the  gravity  model 
is  applied  at  two  different  resolutions  to  determine  whether  or  not  the  number  of  tesserals 
and  sectorials  is  significant.  Also  the  atmospheric  model  and  gravity  models  are  applied 
both  separately  and  together. 

Additionally  all  propagators  demonstrated  a  sinusoidal  nature  in  the  differentials. 
The  commonality  tends  to  suggest  that  it  is  the  data  set  that  is  the  cause  vise  the 
propagators.  Possibly  the  error  is  a  result  of  the  manner  in  which  the  measurements  were 
taken.  An  accelerometer,  for  example,  can  demonstrate  a  slight  sinusoidal  motion  not 
caused  by  external  accelerations. 


Figure  5.4  displays  the  X  Coordinate  Differential.  This  graph  denotes  the 
difference  between  the  predicted  state  vector  and  actual  state  vector  for  the  same  time 
verses  increasing  propagation  time  steps.  Each  differential  is  the  average  result  of 
propagating  all  3238  state  vectors  the  appropriate  time  step  and  taking  the  difference.  As 
expected,  the  Cowell  propagator  shows  the  smallest  differential  even  at  the  largest  time 
step.  The  Encke  propagator  parallels  the  Cowell’s  performance  with  error  most  likely 
due  to  its  assumption  of  perfect  elliptical  motion.  The  Clohessy-Wiltshire  shows  the 
worst  performance.  The  atmospheric  model  produced  no  apparent  effect  which  really  is 
not  that  surprising  considering  the  low  density  of  the  atmosphere  coupled  with  the  short 
time  steps  consider.  The  fact  that  the  addition  of  the  gravitational  model  resulted  in  a 
worse  solution  not  a  better  lends  credence  to  the  possibility  that  the  data  set  contains 
processed  state  vectors  not  raw  state  vectors. 

In  Figure  5.5  the  Y  Coordinate  Differential  is  shown.  All  propagators  show 
almost  identical  error,  including  the  end  point  of  the  unperturbed  Clohessy-Wiltshire. 
This  may  appear  unusual  until  the  magnitude  of  the  Delta’s  is  examined,  which  at  their 
worst  never  exceed  0.3  meters.  The  small  magnitude  of  this  error  is  because  the  motion 
of  both  the  Orbiter  and  Target  spacecraft  are  in  very  nearly  the  same  orbital  plane. 
Considering  the  small  value  of  the  input  component  the  error  is  proportionately  small. 

The  Z  Coordinate  Differential,  shown  in  Figure  5.6,  shows  predictable  results. 
Again  the  Cowell  and  Encke  propagators  are  essentially  error  free.  The  addition  of 
gravitational  and  atmospheric  drag  perturbations  in  the  Clohessy-Wiltshire  produce  a 
negligible  effect.  The  reason  for  the  negligible  difference  when  the  perturbing  forces  are 
added  is  most  likely  due  to  the  magnitude  of  the  velocity  in  the  z-direction.  Since 
atmospheric  drag  in  a  specific  direction  is  a  function  of  the  velocity  of  the  spacecraft  in 
that  direction,  and  because  the  x-axis  of  the  coordinate  system  used  is  defined  by  the 
spacecraft’s  velocity  vector,  drag  can  be  expected  to  have  a  greater  effect  in  the  x- 
direction  than  in  the  z. 
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Figure  5.4  X  Coordinate  Differential 
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Figure  5.5  Y  Coordinate  Differential 
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Figure  5.6  Z  Coordinate  Differential 
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Over  the  given  time  interval,  the  velocity  components  never  exceed  7.5  km/s 
while  the  velocity  differentials  are  all  less  than  2.5  m/s.  Figures  5.7,  5.8  and  5.9  show  the 
Velocity  Differentials  produced.  In  all  three  cases  the  Cowell  propagator  displays  the 
least  error.  The  Encke  propagator  is  generally  the  second  most  accurate,  except  in  the 
case  of  the  Y  Component  where  the  Encke  incurs  it’s  greatest  relative  error.  The  source 
of  this  error  can  be  traced  back  to  the  initial  assumption  made  in  deriving  the  “f  and  g” 
functions,  which  is  that  all  motion  is  assumed  to  occur  in  the  orbital  plane,  also  known  as 
the  x  and  z  plane.  Therefore  the  variation  in  the  component  is  a  measure  of  the  changes 
in  out  of  plane  relative  velocity  of  the  Orbiter  and  Target  occurring  over  the  time  interval. 
This  effect  probably  would  be  seen  in  the  Y  Position  Component  if  the  time  interval  had 
been  greater. 

E.  COMPUTATION  TIME 

The  computation  times  required  for  the  each  propagator  can  be  seen  in  Table  5.1. 
They  were  determined  by  reading  the  computer’s,  biostimer  at  program  start  and  finish 
and  dividing  the  difference  by  the  number  of  vectors  propagated.  As  expected,  the 
Clohessy-Wiltshire  proved  the  quickest  with  Encke  next  and  the  Cowell  taking  the 
longest  time.  Adding  the  perturbation  models  significantly  increase  the  computation  time 
of  the  Clohessy-Wiltshire,  with  little,  if  any,  increase  in  accuracy.  In  fact,  Encke 
propagator  without  any  perturbations  proved  to  be  more  accurate  with  the  second  shortest 
computation  time. 
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Figure  5.7  X  Velocity  Differential 


39 


Time  Step  (s) 


Delta  (km/s) 


o 

o 

o 

b 

O 

o 

b 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

o 

- i. 

o 

N> 

cn 

IO 

cn 

I 


Ht 


he 

e 

O 

< 

£L 

< 

£L 

o 

3 

c r 

n 

W 

o 

3 

i 

o 

?r 

< 

fO 

JL 

o> 

i— 

t 1  + 

£  05  93 

Jj-  *i  >-t 


p 

00 

00 

e4 

•1 

>1 

3 

O 

% 

§ 

e 

i 

u> 

o 

* 

< 

;* 

o 

2L 

w 

o 

'■f 

e 

t 

< 

n> 

< 

G> 

W«-‘ 

Figure  5.8  Y  Velocity  Differentia! 
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Figure  5.9  Z  Velocity  Differential 
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Time  Step  (s) 


seconds/vector 

Cowell 

Clohessy- 

Wiltshire 

Unperturbed 

n/a 

0.0167 

0.0140 

Gravitational  1x0 

n/a 

n/a 

0.0501 

Gravitational 

n/a 

n/a 

0.0796 

30x30 

Atmospheric 

n/a 

n/a 

0.0413 

Atmospheric  & 

0.1303 

n/a 

0.0589 

Gravitational  1x0 

Table  5.1  Propagator  Computation  Times 
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VI.  CONCLUSIONS 


The  analysis  performed  by  this  thesis  has  attempted  to  determine  which  of  three 
representative  orbital  propagators  demonstrates  the  best  performance.  This  performance 
determination  is  based  not  only  on  the  accuracy  of  the  propagator,  but  also  on  the 
computation  time  the  propagator  requires.  Coupled  with  both  of  these  criteria  is  the 
question  of  whether  a  propagator  should  apply  the  extensive  models  required  for 
perturbing  accelerations.  Is  the  gain  in  accuracy  worth  the  increased  time  required  to 
compute  the  results. 

A.  THE  RESULTS 

As  expected  the  Cowell  propagator  demonstrated  the  most  accurate  result  while 
the  Clohessy-Wiltshire  required  the  least  computational  effort.  As  feared  the  perturbing 
accelerations  produced  a  less  accurate  solution.  This  should  be  impossible,  therefore  the 
data  set  used  must  have  been  of  processed  data.  This  fact  does  not  completely  discount 
the  project  findings  because  a  relative  performance  between  the  propagators  can  be 
inferred  based  upon  the  unperturbed  results. 

The  greatest  performance  appears  to  be  demonstrated  by  the  Encke  propagator 
which  without  the  addition  of  perturbing  accelerations  is  essentially  the  “f  and  g” 
functions.  The  Encke  propagator  demonstrated  near  Cowell  accuracy  with  near 
Clohessy-Wiltshire  computation  times.  Only  in  propagation  of  the  y  velocity  component 
did  the  Clohessy-Wiltshire  demonstrate  superior  performance  to  that  of  the  Encke 
methodology. 
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B.  ADDITIONAL  WORK 


In  order  to  fully  complete  this  project  a  number  of  things  need  to  be 
accomplished. 

1.  Raw  Data 

The  entire  thesis  should  be  repeated  using  the  raw  state  vectors  from  a  single 
sensing  source.  This  sensing  source  need  not  be  specified,  but  the  data  set  used  must  not 
have  been  previously  subjected  to  a  propagator  in  order  to  smooth  the  data  set  as  was 
probably  done  in  this  case.  The  C++  code  developed  in  this  thesis  is  already  set  up  to 
perform  this  analysis  whenever  an  appropriate  data  set  can  be  found. 

2.  Encke  Propagator 

While  the  Encke  Propagator  worked  extremely  well  without  perturbing 
accelerations,  this  project  still  requires  a  working  Encke  propagator  with  perturbing 
accelerations  in  order  to  determine  the  potential  benefits  of  adding  these.  The  core  work 
and  theory  exists,  but  the  difficulties  in  compiling  and  executing  the  C++  code  must  be 
resolved.  It  is  not  anticipated  that  the  addition  of  perturbing  accelerations  will  result  in 
any  significant  improvement  in  the  algorithms  performance  for  the  short  time  interval 
calculations.  This  is  especially  true  if  one  considers  the  fact  that  addition  of  the 
perturbing  accelerations  significantly  increases  computation  time. 

3.  Repeated  Iterations 

The  propagators  should  also  be  modified  to  check  for  the  accuracy  of  the 
algorithms  over  longer  propagation  time  in  order  to  define  the  propagator  limits.  This 
could  probably  be  done  by  breaking  the  propagation  time  into  a  series  of  sufficiently 
small  time  steps  to  maximize  performance.  This  is  identical  to  the  method  used  in  the 
Cowell  propagator.  In  this  method  the  state  vector  would  be  propagated  over  a  series  of 
time  steps  with  the  resulting  state  vector  being  used  in  successive  iterations.  The  stability 
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of  the  algorithms  as  well  as  the  additive  effects  of  error  should  be  critically  examined. 
Should  an  operational  algorithm  be  developed  using  one  of  the  propagators  reviewed  in 
this  thesis,  such  an  analysis  should  be  done  in  the  selection  process. 

4.  Different  Orbits 

This  and  the  above  projects  should  also  be  repeated  for  different  orbits  to 
determine  which  performs  best  under  those  circumstances.  More  detailed  knowledge  is 
needed  of  the  limits  at  which  the  perturbing  forces  are  required.  In  this  thesis,  some  of 
the  possible  perturbations  were  discarded  purely  based  on  logical  assumptions.  These 
could  be  examined.  Additionally  the  perturbations  used  probably  have  limits  which 
should  be  examined. 
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